Das Euler-Verfahren

TODO

Lernziele

Lernziele

• Sie können eine Bewegung als Anfangswertproblem formulieren.

• Sie können das Euler-Verfahren zur numerischen Lösung eines Anfangwertproblems in Python implementieren.

• Sie wissen, welchen Einfluss die Schrittweite auf den Fehler des Euler-Verfahrens hat.

Geschwindigkeit aus der Beschleunigung rekonstruieren

Wir führen das folgenden Experiment aus. Eine Person setzt sich auf ein Bobby Car und fährt eine kurze Zeit, so schnell sie kann.

Fig. 1 Symbolbild der Bobby-Car-Fahrt (Quelle: Freepik, Lizenz: Freepik Premium Lizenz)

Zu Beginn der Bobby-Car-Fahrt ist die Geschwindigkeit Null. Wir bezeichnen die Geschwindigkeit mit \(v\) [m/s] und die Zeit mit \(t\) [s]. Es gilt also für den Start \(t_0 = 0 \,\mathrm{s}\)

\[v(0) = 0 \frac{\mathrm{m}}{\mathrm{s}}.\]

Die Beschleunigung messen wir mit Hilfe eines Smartphones und der App Phyphox. Beispielhaft erhalten wir die folgenden Beschleunigungsmessungen, die wir mit dem Python-Modul Pandas importieren.

import pandas as pd

daten = pd.read_excel('data/Beschleunigung_V01.xls')
daten.info()

<class 'pandas.core.frame.DataFrame'>
RangeIndex: 984 entries, 0 to 983
Data columns (total 5 columns):
 #   Column                         Non-Null Count  Dtype  
---  ------                         --------------  -----  
 0   Time (s)                       984 non-null    float64
 1   Linear Acceleration x (m/s^2)  984 non-null    float64
 2   Linear Acceleration y (m/s^2)  984 non-null    float64
 3   Linear Acceleration z (m/s^2)  984 non-null    float64
 4   Absolute acceleration (m/s^2)  984 non-null    float64
dtypes: float64(5)
memory usage: 38.6 KB

Wir verwenden das Python-Modul Plotly Express, um die Beschleunigungsdaten in x-Richtung abhängig von der Zeit zu visualisieren. Die x-Richtung wird gewählt, weil das Smartphone beim Experiment entsprechend montiert war (siehe Koordinatensystem von Phyphox).

import plotly.express as px

fig = px.scatter(daten, x = 'Time (s)', y = 'Linear Acceleration x (m/s^2)')
fig.show()

Damit wir einfacher auf die Messdaten zugreifen können, speichern wir die Zeitmessungen Time (s) in der Variablen t und die Beschleunigung Linear Acceleration x (m/s^2) in der Variablen a.

t = daten['Time (s)']
a = daten['Linear Acceleration x (m/s^2)']

Nun wollen wir numerisch die Geschwindigkeiten aus den Beschleunigungsdaten rekonstruieren. Generell führt eine konstante Beschleunigung \(a\), die eine Zeitdauer \(\Delta t\) lang wirkt, zu einer Geschwindigkeitsänderung \(\Delta v\):

\[\Delta v = a \cdot \Delta t.\]

Da wir nur diskrete Messwerte vorliegen haben, nehmen wir an, dass die Beschleunigung \(a_0\) konstant im Zeitintervall \([t_0, t_1]\) gewirkt hat und können damit nun die neue Geschwindigkeit \(v_1\) berechnen. In mathematischen Formeln erhalten wir

\[v_1 - v_0 = a_0 \cdot (t_1 - t_0)\quad \Rightarrow \quad v_1 = v_0 + a_0 \cdot (t_1 - t_0).\]

In Python ergibt sich mit dem Anfangswert v0 = 0

# Anfangswert
t0 = 0
v0 = 0

# neue Geschwindigkeit
v1 = v0 + a[0] * (t[1] - t0)
print(f'v1 = {v1}')

v1 = 0.0

Jetzt, wo wir \(v_1\) berechnet haben, können wir erneut annehmen, dass \(a_1\) konstant im Zeitintervall \([t_1, t_2]\) gewirkt hat und erhalten:

\[v_2 - v_1 = a_1 \cdot (t_2 - t_1)\quad \Rightarrow \quad v_2 = v_1 + a_1 \cdot (t_2 - t_1).\]

In Python ergibt sich

v2 = v1 + a[1] * (t[2] - t[1])
print(f'v2 = {v2}')

v2 = -0.00035676969787421024

Diese Prozedur können wir immer weiter fortsetzen, bis wir die Geschwindigkeit auf dem gesamten gemessenen Zeitintervall rekonstruiert haben. Es ist einfacher, Python die Arbeit machen zu lassen. Daher bestimmen wir zunächst die Anzahl der Messwerte \(N\) und initialisieren anschließend ein NumPy-Array der Länge \(N\), das die Geschwindigkeitswerte aufnehmen soll.

import numpy as np

N = len(t)
v = np.zeros(N)

Die obige Berechnung der Geschwindigkeitsänderung

\[\Delta v = a \cdot \Delta t\]

lässt sich für beliebige Zeitintervalle \([t_i, t_{i+1}]\) formulieren als

\[v_{i+1} - v_{i} = a_{i} \cdot (t_{i+1} - t_{i})\quad \Rightarrow \quad v_{i+1} = v_{i} + a_{i} \cdot (t_{i+1} - t_{i}).\]

Wird eine for-Schleife verwendet, berechnet der folgende Python-Code nacheinander die Geschwindigkeitswerte aus den diskreten Beschleunigungsmessungen:

for i in range(1, N-1):
    v[i+1] = v[i] + a[i] * (t[i+1] - t[i])

Nun können wir die berechneten Geschwindigkeiten visualisieren.

fig = px.line(x = t, y = v,
  labels = {'x': 'Time (s)', 'y': 'Velocity in x (m/s)'})
fig.show()

Anfangswertprobleme allgemein mit dem Euler-Verfahren lösen

In dem obigen Beispiel haben wir aus der Messung einer Beschleunigung die Geschwindigkeit rekonstruiert. Allgemein betrachtet kann der Zusammenhang zwischen Geschwindigkeit und Beschleunigung über die Differentialgleichung

\[\dot{v}(t) = a(t)\]

beschrieben werden. Es liegt eine Differentialgleichung 1. Ordnung vor. Zusammen mit der Forderung, dass zu einem reellen Startwert \(t_0\) die Geschwindigkeit gleich der Anfangsgeschwindigkeit \(v_0\) ist, also

\[v(t_0) = v_0\]

gilt, stellt das Bobby-Car-Problem ein Anfangswertproblem dar. Wir können die oben beschriebene Prozedur auf beliebige Anfangswertprobleme übertragen.

Wenn das allgemeine Anfangswertproblem

\[\dot{x}(t) = f(t, x(t)), \; x(t_0) = x_0\]

lautet, dann beginnen wir beim Anfangswert \(t_0\) und nutzen die Schrittweite \(h > 0\), um Zeitpunkte

\[t_{i} = t_0 + i\cdot h, \quad i = 0, 1, 2, \ldots\]

zu erzeugen. Diese Punkte werden als Gitter bezeichnet. Dann wird die Lösungsfunktion \(x\) durch die Werte

\[x_{i+1} = x_{i} + f(t_{i}, x_{i}) \cdot h, \quad i = 0, 1, 2, \ldots\]

angenähert. Dabei haben wir die Abkürzung \(x_{i}\) für \(x_{i} = x(t_{i})\) verwendet. Wir hangeln uns also immer ein Stückchen “h” weiter. Diese Prozedur wird Euler-Verfahren genannt.

Begründen lässt sich das Euler-Verfahren mit der Taylor-Formel bzw. der Tangentengleichung. Das Taylorpolynom vom Grad 1 an einem Entwicklungspunkt \(t_{*}\) lautet

\[f(t) = f(t_{*}) + f'(t_{*})\cdot (t-t_{*}).\]

Übertragen auf unser Anfangswertproblem gilt für den Entwicklungspunkt \(t_i\) die Approximation bzw. Taylorploynom

\[x(t) = x(t_{i}) + \dot{x}(t_{i}) \cdot (t-t_{i}).\]

Jetzt können wir für die erste Ableitung \(\dot{x}\) die rechte Seite der Differentialgleichung einsetzen und erhalten

\[x(t) = x(t_{i}) + f(t_{i}, x(t_i)) \cdot (t-t_{i}).\]

Speziell für den nächsten Gitterpunkt \(t_{i+1}\) gilt dann

\[x(t_{i+1}) = x(t_{i}) + f(t_{i}, x(t_i)) \cdot (t_{i+1}-t_{i}),\]

was mit den Abkürzungen \(x_i=x(t_i)\) und \(x_{i+1} = x(t_{i+1})\) vereinfacht geschrieben werden kann als

\[x_{i+1} = x_{i} + f(t_{i}, x_i) \cdot h.\]

Einfluss der Schrittweite

Wie gut die numerische Lösung eine exakte Lösung eines Anfangwertproblems annähert, wird auch durch die Schrittweite bestimmt. Um den Zusammenhang zwischen der Schrittweite \(h\) und den Fehler des Euler-Verfahrens zu verstehen, betrachten wir das folgende Beispiel.

Ein Bobby-Car fährt mit einer Anfangsgeschwindigkeit von \(v_0 = 5 \mathrm{m}/\mathrm{s}\). Der Fahrer oder die Fahrerin lässt das Bobby-Car ausrollen. Dann gilt für das Abbremsen des Bobby-Cars die folgende Differentialgleichung

\[\dot{v}(t) = - k \cdot v(t)^2,\]

wobei wir annehmen, dass der Abbremsvorgang zum Zeitpunkt \(t = 0\,\mathrm{s}\) beginnt. Die Konstante \(k\) ergibt sich aus dem Luftwiderstand und hat in unserem Fall einen Wert von

\[k = 0.003 \frac{1}{\mathrm{m}}.\]

Als nächstes variieren wir die Schrittweite \(h>0\).

Mini-Übung

Öffnen Sie im Browser die interaktive Demonstration des Euler-Verfahrens

</Users/gramschs/Bewerbung/book_thma/doc/_static/assets/interactive_euler_demo.html>

für die Abbremsung des Bobby-Cars. Variieren Sie die Schrittweite h und beobachten Sie, wie gut oder schlecht die Annäherung an die exakte Lösung ist. Legen Sie eine Tabelle an:

Schrittweite h	Fehler
1	?
5	?
…	???

Stellen Sie eine Vermutung an. Wie hängt der Fehler von der Schrittweite ab?